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Abstract — The Green’s function for the transport of ions of 
high charge and energy is utilized with a nuclear fragmenta- 
tion database to evaluate dose, dose equivalent, and RBE for 
C3H10T1/2 cell survival and neoplastic transformation as a 
function of depth in soft tissue. Such evaluations are useful to 
estimates of biological risk for high altitude aircraft, space 
operations, accelerator operations, and biomedical applica- 
tions. 
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INTRODUCTION 

The biological response of living tissues depends (in 
part) on the temporal and spatial fluctuations of the 
energy deposits of the ions within the tissue system. Such 
fluctuations depend not only on the specific environment 
(or particle fluence) to which the person is exposed but 
also on how that environment (or particle fluence) is 
modified by interaction with the human body in reaching 
specific tissues. Only by knowledge of the specific 
radiation types and their physical properties at the tissue 
site can a basis for estimating risk be found. Even if the 
environment to which the individual is exposed is known 
precisely, the energy deposits within specific tissues deep 
in the body are largely known through theoretical esti- 
mates and are therefore limited by the uncertainty in the 
calculational models. Although methods for exposure 
evaluations in high energy nucleonic fields are relatively 
well developed (Alsmiller et al. 1970; Nabelssi and 
Hertel 1993), there are several practical applications in 
which exposures with high charge and energy (HZE) 
ions are of concern. Among these are the natural cosmic 
ray environment consisting of relativistic nuclei of all 
elements for which HZE ions contribute 80% of the free 
space dose equivalent (Wilson et al. 1991). Even for 
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modest depths of shield penetration (30 g cm 2 ), the 
HZE dose equivalent is on the order of 45% of the total 
and is of considerable importance to protection of astro- 
nauts and possible future high speed aircraft (^ Mach 2.5 
or greater) on high latitude routes (for example, New 
York-London, London-Anchorage, Frankfort-Tokyo). 
Computational models are required for treatment planning 
in which cell killing for HZE tumor therapy is the critical 
parameter. With the further development of HZE accelera- 
tors, there is a growing concern among policy setting bodies 
on evaluation of exposures for the scattering of HZE 
particles from the beamline (L. Brackenbush 1994 § ) for 
which the present calculations were undertaken. Clearly, an 
accurate conversion of the environment to estimates of 
exposure fields at specific tissue sites is an important issue 
in HZE radiation protection problems. 

In our previous paper on the closed form solution for 
the HZE Green’s function, we reviewed the computa- 
tional procedures for space and laboratory HZE transport 
of various groups (Wilson et al. 1990). It is noted in that 
review that numerical solution methods for the Boltz- 
mann equation are best suited to space radiations where 
energy spectra are smooth over large energy intervals and 
less suited to the simulation of laboratory beams which 
exhibit large spectral variation over a very limited energy 
domain and a large energy derivative resulting in large 
truncation errors (Wilson et al. 1990, 1991). In the 
present paper, we apply these analytical methods (Wilson 
et al. 1990) to evaluate dosimetric quantities as a function 
of penetration depth in tissue equivalent materials to 
estimate specific organ exposures. 

TRANSPORT EQUATION 

The Green’s function is introduced as a solution of 
the following transport equation: 

dx ~ df VO + (T I G jm{x, E, E') = 'E<T jk G k Jx, E, E‘), 

-* k 

( 1 ) 

where G jm {x, E , E f ) is the flux of ion type j at depth x 
with energy E (in MeV/amu) resulting from a unit flux of 

§ Personal communication, L Brackenbush, Pacific Northwest 
Laboratory, POB999, Richland, WA 99352, 1994. 
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type m ions of energy E' at the boundary, Sj(E) is the 
change in E per unit distance, a } is the total macroscopic 
reaction cross section, and Uj k the macroscopic cross 
section for collision of ion type k to produce an ion of 
type j. The boundary condition is given by 

G jm { 0, E, E r ) = 8 jm 8(E - £'), (2) 

where 8- is Kronecker’s delta and 8(E~E') is Dirac’s 
delta. The dose and dose equivalent for monoenergetic 
ions are given by 

Djx, £') = 2 I AjSj(E)GjJx, E, £') dE; and (3) 

J 

H m (x, E') = 2 AjSj(E)Q(Lj)G jm (x, E, E') dE. 

j J 0 

(4) 

If G Jm (x, E , E f ) is known as a closed form solution 
(Wilson et al. 1990), then eqns (3) and (4) can be 
evaluated by simple integration techniques and the asso- 
ciated errors in numerically solving eqn (1) are avoided 
(Wilson et al. 1991). 

The above equations can be simplified (Wilson 
1977) by transforming the energy into the residual range 
of the ions as 


\ E de 

Sj(eY 


(5) 


and defining new field functions as 

r p r’ m ) = Sj(E)G Jm (x, £, E'). (6) 

As a result, eqn (1) becomes 

d a 


dx dr, 


+ CTj )%Jx, r p O = 2 ~<T jk %Jx, r k9 rj, 


v k 


with the boundary condition 

%m(0, /> r'J = 8 jm 8(r, ~ r' m ). 

Note that Vj is the range scale factor such that vp: 
i \„r m and can be expressed as 

Vj = Zj/A r 


(7) 

( 8 ) 

(9) 


The solution to eqn (1) can be written as a perturbation 
series as 


fp r m ) S Gj m (x 9 r p r rn ). 


( 10 ) 


where 

%!(x, r p r'J = exp( - ajx)8 jm 8(x + rj - rj, (11) 


and the first collision term is 

, i Vj<J pn 

%n {x, rj, rj = z - w j exp( - 1 TjXj - r v„), (12) 

i *7 1 

where <8$ (x 9 r p r ni ) takes on non-zero values for 

v ir m -x)<rj< -~r m -x 9 (13) 

and 

Xj = (v„x + v/j - v m r; n )/(v m - Vj); and (14) 

= ~ Vj X ~ V j rj)/{v m ~ Vj). (15) 

Note that in eqns (12) — (15), Xj + x„ = x for all rj and r’ m . 
The significance of x m is that it is the distance ion m 
traveled from the boundary to the collision site at which 
the ion j was produced and must now travel distance jt . 
before reaching x. Higher order terms in the perturbation 
series are given in the Appendix. 

To this point of development, the perturbative 
Green’s function (eqn 10) can be approximated by 
summing the terms given by eqns (1 1), (A7), and (A9). 
The problem in explicitly summing eqn (10) is that the 
higher order terms of eqn (A9) contain thousands of 
terms, each of which are composed of sums of tens of 
thousands of terms for the second order correction (/ = 2) 
and millions of terms for the third order correction. Since 
the terms given by eqns (A7) and (A9) are dependent 
only on spatial coordinates, we seek an alternative means 
of summing the series (10) implicitly. 

NONPERTURBATIVE GREEN’S FUNCTION 

We now introduce nonperturbation terms for the 
summation in eqn (10). First we recall (Wilson et al. 
1990) that the g-function of n arguments is generated by 
the perturbation solution of the transport equation ne- 
glecting ionization energy loss given by 

{ dx + O'; j gjm(x) = 2 (Tjkgkmix), ( 1 6) 

subject to the boundary condition 

gjM = 8 jm (17) 

The solution is 

gjjx) = 8j m g(m) + a Jm g(j, m) + 2 o- jk (r hn g(j, k,m) + ... 

(18) 

A term by term comparison of the series (18) with the 
series (10) shows that in the spectral average approxima- 
tion (see Appendix) we may write 

rp r, n ) ' exp( cTjx)8j m 8{x T k. r m ) 

( 19 ) 

v,[gjm(x) - exp( - <jjx)8jJ 

x( v m - Vj) 
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The advantage of eqn (19) is that g jm {x) satisfies the 
convolution product relation (Wilson et al. 1993a) as 

gjmix) = S gjk( X ~ y^kmiy) (20) 

k 

for any positive values of x and y. A solution for small x 
is easily obtained from the first few terms of eqn (18). 
Eqn (20) may then be used to propagate the function 
gj nt (x) over the entire solution domain taking arbitrarily 
large steps. Although the nonperturbation Green’s func- 
tion in the form of eqns (19) and (20) provides a rapid 
computational method, the spectral terms are replaced by 
averages over the spectral domains of higher order terms. 
These energy averaged spectra can be corrected with the 
spectral distributions of perturbation theory. The first 
collision Green’s function is given as 


ty„(x, /> r’ m ) = r | exp( - 1 TjXj - <r„ r xj. (2 1 ) 

We rewrite in terms of its average value (given by 
eqn A7) as 


^ Vjr r . ) — 

yJ jm r p ' mr 

where (x) 


g)m( x ) 


+ bjm(x)(rj - rj), (22) 


(v„, - Vj) X 
— <T /fn g(j\ m) is given by eqn (A7) and 


r ji + r ju 


is the midpoint r y between its limits given by eqn (A8). 
The b jm term of eqn (22) has the property that 


I b jm {x)(r ~ Fj) dr = 0 (23) 

^ r„ 

to insure that the first term of eqn (22) is indeed the 
average spectrum as required. The spectral slope param- 
eter is 


i^Jexp) — ayx) - exp( — cr„,.r)] 
b jm (x) = — 7 ^ j ■ (24) 

' X{V m - Vj) j V m - Vj\ 


Note that the spectral slope results from the difference in 
mean free paths of the projectile and secondary ion. The 
ion j at .v produced near x m = 0 is less attenuated than 
those produced locally (x m = x). Since those produced 
locally are generally of lower energy than those produced 
at x = 0 we obtain the spectral slope. 


DEPTH-DOSE RELATIONS 

We consider the evaluation of dose equivalent as 
given by eqn (4) using the present formalism. The quality 
factor chosen is the latest recommended values of the 
ICRP (1990) given by 

[ 1 L<10keV//jun 

Q(L) = 0.32L — 2.2 10 < L < 100 keV//xm (25) 

[ 300/ yjL L > 1 00 keV/jxm 


It is clear from the above formalism that 
HJx, E’) = Q[L m (E 0 )]Aj m (E a )e-™ 

VjigjJx) - S >m exp(-cr„, jc)] 


+ 'ZAf 


(k, - Vj)x 


fij- 

X Q[Lj(E)] dE (26) 

k' ; ! 

f 

+ 2 AjbjJx) 0[L/£)][/?/£) - #j] dE 

j b Jt 

where E ju and E j{ are obtained from the limits in eqn (13). 
D m (x, E') is similarly obtained from eqn (26) by setting 

G(t)=> l. 

The ion flux within a water column is represented by 
the eighty most important isotopes lighter than 59 Ni, The 
Green’s function contains 6,400 terms as related by 
3,160 fragmentation cross sections. The ion flux within 
the column, due to a monoenergetic incident ion flux, is 
evaluated using the NUCFRG database (Wilson et al. 
1987; Townsend et al. 1993; Wilson et al. 1994) for the 
oxygen constituents and the Silberberg-Tsao (Silberberg 
et al.1976) model for the hydrogen constituents. The 
isobaric nuclear fragmentation cross sections are shown 
in Fig. 1. The nuclear absorption cross sections are 
shown in Fig. 2 for selected ions. In the following 
calculations, the nuclear cross sections were evaluated at 
the ion beam energy at the entrance to the water column. 
The results of the computations using eqn (26) with Q 
1 for D m {) c, E’ ) are shown in Fig. 3 for 20 Ne beams and 
40 Ar beams in a water column with the experimental 
measurements of Schimmerling et al. (1989) and Lyman 



1 

D 


(c) 600 MeV/amii. 


(d) 2400 MeV/nmii. 


Fig. 1. Isobaric fragmentation cross sections projectiles of 6 Li to 
sy Ni for four incident energies in water targets. 
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Energy (MeV/nucleon) 

Fig. 2. Nuclear absorption cross sections of selected projectiles in 
water. 


et al. (1977). Reasonable agreement is obtained in each 
case although there are important differences. These 
differences are on the order of 20% or less and are 
consistent with the comparison of LET spectral measure- 
ments of specific fragments (Shavers et al. 1993) mea- 
sured in the 20 Ne experiments. Improvement must come 
through more accurate nuclear models, addition of en- 
ergy dependence in the nuclear cross sections, and 
contributions from target breakup. 

The average quality factors H m ( jc, E')/D m (x , E f ) for 
fixed range for beams of 20 Ne, 40 Ar, and ^ 6 Fe are shown 
in Fig. 4. The quality factor is relatively constant over the 
range of the primary beam except in the Bragg peak 
region. Beyond the Bragg peak the average quality factor 
for "°Ne beams drops suddenly as expected and declines 
slowly to larger depths since only the lighter ions 
penetrate to large depths. The depression of the average 
quality factor at the Bragg peak of the 40 Ar and 56 Fe 
beams is due to the decline of the primary ion quality 
factor ( Q «=* 300/VZ) and the rapid rise in the primary 
ion dose at the Bragg peak. 

There is experimental evidence that relative biolog- 
ical effectiveness depends not only on the rate at which 
an ion gives up energy to the tissue medium (numerically 
equal to LET) but depends as well on the lateral extent of 
the energy deposit (track structure effects). For this 
reason we consider a track structure cellular repair model 
for survival and neoplastic transformation (Wilson et al. 
1991) as applied to the C3H10T1/2 mouse embryo cell 
culture for which there are extensive data using HZE ion 
beams (Yang et al. 1985, 1989). The cellular repair 
model is derived from the following assumptions. Injury 
from gamma rays follows Poisson statistics with a 
characteristic dose D (> and enzymatic repair is less 
efficient with increasing number of hits. Such a model is 
consistent with the concept of a high efficiency and fast 
but saturable repair enzyme pool in competition with a 
slower less efficient repair process (Tubiana et al. 1990). 
Injury within a particle track is mediated by the second- 
ary electrons (delta rays) and is related to the gamma-ray 




Depth in water (cm) 

(b) Argon 

Fig. 3. Bragg curves for (a) Ne and (b) Ar ion beam according to 
present calculations and experiments. 

injury through An "inactivation 1 cross section asso- 
ciated with the core of the ion track is approximately 
given by 


o- = cr fl (l - e - z ' 2 /K(} 2 r, (27) 


where k is a nondimensional size parameter associated 
with the sensitive structure within the cell, Z is the ion 
effective charge, and /3 is the ion velocity in units of the 
velocity of light (Katz et al. 1971). In the track periphery, 
the injury level is below saturation and multihit repair 
kinetics as in the case of gamma-ray exposures are 
assumed and the probability of direct inactivation is 
given as p = c r/cr () where cr () is the Katz ‘saturation 1 cross 
section. The probability of injury with multihit repair 
kinetics is (1 - P). In the low dose rate limit (Wilson et 
al. 1993b), the nonsurviving fraction of exposed cells is 


>h„(t) a m] (1 - P)D a 
= — 6 1 ' 3 - - - ,D 

n„ a, D„ L 


( 28 ) 
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(b) Argon (620A MeV) 



where D is the accumulated dose, L is the linear energy 
transfer, a ml is the misrepair rate for once hit cells, a x is 
the total rate at which the enzyme (Z) forms a repair 
complex (C') and completes the repair through the 
reaction kinetics given by 

Z Z + DNA 
Z + DNA, -» C N z + DNAm 

where subscript / denotes injured DNA and subscript m 
denotes misrepair. Note that (1 - a m[ ia x ) and a^/a, are 
the branching ratios of the reaction to perfectly repaired 
and misrepaired states. In eqn (28), the 6 1/3 (1 —P)D/D {) is 
the fraction of cells with injury in the track periphery and 
ct m JcL is the probability the cell is permanently harmed. 
The fraction of cells exposed within the ion core with no 
chance of recovery is aD/L. A similar result holds for the 
fraction of transformed cells in the low dose rate limit 
with appropriate parameters for transformation. The 
repair rates and efficiencies (which depend on the cell 
cycle status) are found from the experiments of Yang et 
al. (1985, 1989) and the low dose rate limit parameters 
for resting G 0 phase and exponential growth phase cells 
are given in Table 1. Unlike conventional dosimetric 
analysis wherein radiation quality is represented by LET 
dependent factors, the repair kinetics model is driven by 
track-structure dependent injury coefficients (eqn 27). 
The RBE W is found by finding the gamma-ray dose and 
ion dose which leads to the same fraction of changed 
cells (note P = 0 for gamma rays) as 


RBE W = 1 - P + 


a | o-D 0 
®-m i 3 


(29) 


The subscript m on RBE denotes the maximum value 
obtained in the low dose rate limit. The RBE, N values 
derived from a cell kinetics model (Wilson et al. 1993b) 
for C3H10T1/2 cells exposed and repaired in stationary 
phase at the low dose rate limit are substituted for Q(L) 
in eqn (27) for ion beam exposures. The RBE m for cell 
survival (S) and cell neoplastic transformation (T) are 
shown in comparison to the average quality factors. The 
quality factor is intended to represent stochastic pro- 
cesses (mainly cancer induction) and should best corre- 
late with the RBE for neoplastic transformed cells (a 
pre-cancerous cell). The RBE for cell survival is more 
indicative of deterministic effects for which large RBEs 


Table I. Cellular kinetic ratios and parameters. 

Cellular kinetic ratio ct n ^ot x - \ 

G 0 phase Exponential phase 

Survival M1.03 0.30 

Transformation .002 .01 


Katz C3H10T1/2 cell parameters 

cr 0 , cm 2 k m D u , Gy 

5X10" 7 750 3 2.8 

7X10 " 475 3 116 


Fig. 4. Average quality factor and RBE,„ values forC3H101/2 cell 
survival (S) and transformation (T) for ions of 20 cm range in 
water; (a) 2() Ne, (b) 40 Ar, (c) 56 Fe. 


Survival 

Transformation 
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have been measured for mouse embryo hemopoiesis 
(Jiang et al. 1994) irradiated with high LET alpha 
particles. Many of the qualitative features of the RBE m 
are represented by the average quality factor but impor- 
tant track structure dependent differences can be ob- 
served. The quality factor increases from beam entrance 
to the position at which the primary ion reaches 100 
keV/juim and declines to greater depth as seen for neon 
and argon ions. The iron ions enter with LET >100 
keV/jLim and the quality factor declines throughout the 
medium. The quality factor approaches the stochastic 
RBE(T) for the lighter more penetrating fragments 
downstream from the Bragg peak. Near beam entrance, 
the width of the primary track is maximum and the 
RBE(T) lies substantially below the LET dependent 
quality factor. It is accidental that RBE(S) coincides with 
quality factor at beam entrance for neon and argon. The 
diverse behavior of RBE(S) and RBE(T) is a result of the 
specific track structure parameters associated with sur- 
vival and transformation (Table 1) and are consistent 
with the data of Yang et al. (Wilson et al. 1993b). The 
quality factor is most appropriate for representing, the 
cell transformation RBE W for lighter ions such as ~°Ne 
and appears to be conservative for heavier ions. A quality 
factor for deterministic effects is not defined but would 
be more indicative of RBE(S). Note that the present 
result neglects the secondaries produced from tissue 
constituents (target fragments). 

CONCLUSIONS 

At this point of development of HZE transport 
theory, the fluence spectra of the multiple charged ions 
can be reasonably calculated in one dimension assuming 
constant nuclear cross sections. Clearly, future efforts 
should concentrate on deriving important corrections 
associated with realistic cross section variation with 
energy. Even with this improvement several issues re- 
main in comparing the results to the experimental data. 
The target breakup will make contributions at all depths 
but may be a source of error near the Bragg region and 
beyond and especially for the lighter ion beams. Beyond 
this, significant improvements are required in the nuclear 
models used in generating the nuclear database. 
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APPENDIX 


The second collision term is similarly 


V \ f- p* ) = 

V -* •> ' n 1 m) 


2 


Vj(r jk (r km 

V m - V k \ 


(Al) 


X exp( — cr,.*-, - (T k x k - a m xj d.\j, 


where 


Xj + x k + x m = x (A2) 

+ *vv„, = v„/ m - v k r h (A3) 

and the integral limits in eqn (Al) depend on v m9 v k , v f 
and the physical requirements that the distances x j9 x k , 
and x m all be positive. We have used the rather weak 
dependence of the exponential argument of eqns (1 2) and 
(Al) on energy to simplify spectral quantities in the past 
(Wilson et al. 1984). Higher order terms are similarly 
derived (Wilson et al. 1990). Whenever r j9 r’ m > x it can 
be shown that eqn (12) has the property that 

f v 

r p r’J dr, = <r jm g(j 9 m) 9 (A4) 


* * * *jn>jn+\) 

(A6) 

• • * ijn-uj/t) g(j 1 5 7*2, '•••>]> i 1 ijn+i) 


The first collision Green’s function term may be approx- 
imated by the spectral averaged value as 


* ■ ( r r r f ) 

' ]i r in/ 


%m(x 9 r j9 rj drj 


(A7) 


_ Vj<T jm g(j, m) 

x{ v m - Vj) 

and is functionally dependent only on depth of penetra- 
tion x (Wilson and Badavi 1992). Its contribution applies 
only over the energy range 


(r, + x) - C < 1 r } + x ( A8) 

"m 

for v m > v-. If v f > v m9 as can happen in neutron removal 
collisions, the negative of eqn (A7) is used and the 
inequalities of eqn (A8) are reversed (Wilson et al. 1990). 
The higher order terms are approximated (Wilson and 
Badavi 1992) by their spectral averages as 


where g(y, m) is the g-function of two arguments intro- 
duced in our prior paper (Wilson et al. 1990) and the 
integral limits are the interval limits defined by eqn (13). 
The g-functions are given by the following relations 

g{j) = cxp(-oyx) (A5) 


and 


<€(*, n, r;„) 


1 

k u k: ki 


v , (T jh (T hh' • <r k k u k 2 , ■ ■ % k, _|. m) 

X( V,„ ~ Vj) 


(A9) 

Note that (x, r r r,',) is purely dependent on .v for / > 
0, and we represent the expression as %[',!( •*)• 




